Fragility of the Schrödinger Cat in thermal environments

We describe the decoherence instability of Schrödinger Cat states in the two-site Bose-Hubbard model with an attractive on-site interaction between particles. For N particles with onsite attractive energy U and hopping amplitude between sites t, Cat states exist for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u\equiv \frac{UN}{2t}<-1$$\end{document}u≡UN2t<-1 at zero temperature. However, they are increasingly unstable to small thermal fluctuations as the Cat itself is increasingly well-defined and its components become well-separated. For any given \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u<-1$$\end{document}u<-1, the decoherence temperature becomes smaller for large N. The loss of off-diagonal coherence peaks in the equilibrium density matrix is dominated by the thermal admixture of the first excited state of the many-body system with its ground state. Particle number fluctuations, described in the grand canonical ensemble also reduce coherence, but to a lesser degree than thermal fluctuations. The full density matrix of the Schrödinger Cat is obtained by exact numerical diagonalization of the many-body Hamiltonian and a narrow regime in the parameter space of the particle number, temperature, and U/t is identified where small Cat states may survive decoherence in a physical environment.


Analytical model
The boson number operator commutes with the Hamiltonian [H, N] = 0 .Thus, energy eigenstates can be chosen to have a definite total number of bosons.We diagonalize the Hamiltonian in equation (1) using the orthogonal basis of Fock states 37,44 : | � = l l |l� , where |l� represents the normalized quantum state |l, N − l� describing l bosons on site a and N − l bosons on site b.Here, l takes values from 0 to N. The coefficients l are yet-to-be- determined wavefunction amplitudes.After substituting this ansatz into the Schrödinger equation corresponding to the Hamiltonian (1), we obtain the eigenvalue equation: where + UN 2 ( N 2 − 1) and E is the energy eigenvalue.This describes a quantum particle hopping in a one-dimensional lattice with N + 1 sites, in a potential V l and for which the hopping amplitudes vary with location, with hopping toward either endpoint of the lattice becoming increasingly difficult.We perform numerical diagonalization of the above Hamiltonian with an orthonormal Fock state basis of N + 1 vectors (N bosons), using the Mathematica software 45,46 .We discuss the exact numerical results later.
A more intuitive physical picture is obtained by mapping the above lattice problem to an effective continuum model.We show, below, that the finite difference equation (2) can be mapped to a Sturm-Liouville differential equation in a continuous variable, y, that ranges from 0 to N. We introduce the definition of first difference �� �l | l+1/2 ≡ � l+1 − � l .The above eigenvalue equation (2) can then be expressed in the form of a discrete Sturm- Liouville problem: E� l = − � �l t l−1/2 ) is an effective lattice potential.We extend the discrete wavefunction amplitude l to a continuous and differentiable function �(y) .According to the Lagrange mean value theorem, d� dy (y = ζ ) = � l+1 − � l for some ζ in the open interval (l, l + 1) .Using the ansatz, ζ = l + 1/2 , the eigenvalue equation in continuum limit becomes: E�(y) = − d dy t y−1/2 d�(y) dy − t y − t y−1 + UN 2 ( N 2 − 1) and �(y) represent the effective potential and the wavefunction, respectively, in the continuum limit, where 0 ≤ y ≤ N .Here, t y = t (N − y)(y + 1).
For convenience, we introduce the translated variable x ≡ y − N/2 , define over the symmetrical interval The eigenvalue equation is expressed in terms of x: where t x = t N 2 N 2 + 1 − x(x + 1) can be interpreted as the reciprocal of the position-dependent effective mass in the continuum model. (1) www.nature.com/scientificreports/

Effective potential in the continuum model
We identify the critical transition by the change of concavity of V eff (x) , and map the second order differential equation (3) to that of a simple harmonic oscillator (SHO).The effective potential V eff (x) can be expressed in Taylor series about x = 0: where The coefficients of all non-vanishing terms in the series are positive, except for the x 2 and x 0 terms.The coefficient of . For large N, the effective potential V eff (x) is simplified to Under the simplifying condition, the coefficient of x 2 changes its sign at UN 2t (≡ u) = −1 .We refer to u = −1 as a critical point.For u ≥ −1 , V eff (x) has a single minimum at x = 0 .For u < −1 , V eff (x) exhibits a double well structure with a pair of local minima, symmetric about the origin.Figure 1a shows the behaviour of the effective potential for two different values of u.The effective potential at u = −0.5 , which is above the critical value, is shown by the blue line.The green line for u = − 1.25 is below the critical value and exhibits two separated minima.
When u < −1 , the double well potential minima occur at x = ±x c (see Fig. 1c).For large N, we set dV eff /dx(x) = 0 to determine In the vicinity of the local minimum at x c , the effective continuum model can be approximated as a simple harmonic oscillator of mass m c and angular frequency ω c ( = 1): To complete the analogy with a SHO, we identify m c using the kinetic energy term in Eq. ( 3) and ω c using the quadratic term of V eff (x) at x = x c in Eq. ( 8): By direct comparison with SHO, we estimate the ground state energy E g ≈ (1/2)ω c relative to the potential minimum V eff (x = x c ): (5) In this analogy, the ground state solution is proportional to exp − (x , where the width Similarly, the left side minimum has a ground state proportional to exp − (x . Given the non-zero tunnelling amplitude between the minima, the ground state wave function for the double well potential can be approximated by a linear combination of � (0) R (x) and � (0) L (x) .The numerical solution of the lattice model (2) for the ground state wave function in the double well for u = − 1.25 is shown by the green line in Fig. 1b.
Well above critical point when u = − 1 , the ground state wave function exhibits a single Gaussian peak at origin.Near x = 0 , the effective continuum model can be approximated as a SHO of mass m 0 and angular frequency ω 0 .Similar to the previous analysis, by Eqs. ( 3) and (8), we identify This approximate ground state wave function is proportional to exp − x 2 /(2σ 2 0 ) , with σ 0 ≈ ( N 2 ) 1/2 (u + 1) −1/4 .The numerical solution of lattice model (2) for the ground state wave function at u = − 0.5 is shown by the blue line in Fig. 1b.
A ground state in which the wave form L consists of distinguishable peaks at x = ±x c corresponds to many-body quantum state of (1) that is a superposition of a significant number of bosons on site a with the same significant number of bosons on site b.In what follows, we define a well-developed Schrödinger Cat state to be one in which the ground state probability density at x = ±x c is at least a factor of ten larger than the prob- ability density at x = 0 .As u → −1 , from below, the two components of the "Cat state" merge and the peaks at x = ±x c becomes less discernible.On the other hand, for large negative u < −1 , the tunnelling amplitude between the two components of the "Cat" becomes negligible and we approach an "extreme Cat" state of the form (|N, 0� + |0, N�)/ √ 2 .While such "extreme Cat" states may appear as tantalizing possibilities at zero temperature, we show below that they are unstable to decoherence for infinitesimally small thermal fluctuation.
So far, we have presented analytical approximations describing a single minimum in V eff (x) .We describe below an analytical approximation to the tunnelling amplitude between a pair of distinct local minima that accurately describes the ground and first excited state wavefunctions, provided we are not too close to the critical point u = −1.

Approximation for ground and first excited states
We apply the WKB approximation to the continuum two-site Bose-Hubbard model to estimate the energy separation dE 10 between the ground and first excited states.For u < − 1 , the ground state and first excited state are even and odd superpositions, respectively, of near-Gaussian wavefunctions centered at the minima of the double well potential.This double well potential intersects the ground state energy at four different classical turning points ( ±x 1 and ±x 2 ), as depicted in Fig. 1c.
The WKB approximation for the ordinary Schrödinger's equation is modified to accommodate a positiondependent effective mass in the Sturm-Liouville problem: Here, we introduce ǫ = 1 to keep track of the perturbation order of the WKB approximation.We begin with an exponential asymptotic approximation 47 : where S n (x) is the nth order term of the phase function.Here, we only keep the phase function up to the first order in ǫ .We substitute Eq. (15) into Eq.( 14) to obtain a sequence of equations which determines S 0 and S 1 : where Q(x) ≡ (V eff (x) − E)/t x−1/2 .These equations are solved to yield: www.nature.com/scientificreports/As a result of the position-dependent effective mass, the modified WKB wavefunction gains an additional amplitude modulation factor compared to the WKB wavefunction for the standard Schrödinger's equation 47 : As shown below the extra factor does not affect the energy quantization conditions.The WKB approximation fails near a classical turning point.Nevertheless, the effective potential can be linearized in the vicinity of a turning point, and the solutions to the linearized differential equation are the Airy functions.A global approximation can be constructed by matching the asymptotic expansions of the Airy functions to the WKB approximation on either side of the turning point: a phase angle over the classically allowed region x 1 < x < x 2 .The standard con- nection formulae apply because the effective mass is regular and positive at the classical turning points.Over the classically forbidden region, 0 ≤ x < x 1 , there is a factor of 2 with the cos η term, because of the different prefactors in the asymptotic expansions of the exponentially decaying Airy function Ai(x) of the first kind and exponentially increasing Airy function Bi(x) the second kind.
The WKB approximation for x < 0 is determined by parity.For the even-parity ground state, � WKB (−x) = � WKB (x) .For the odd-parity first excited state, � WKB (−x) = −� WKB (x) .These imply the following quantization conditions: d� WKB /dx(x = 0) = 0 for even parity states and � WKB (x = 0) = 0 for odd parity states.Hence, by Eq. ( 19), we obtain where κ = a phase angle in the non-classical region.Here, + and − represent even and odd parity states, respectively.The calculations are simplified by observing that t x−1/2 and Q(x) are even functions of x, so that their derivatives vanish at x = 0.
Finally, we estimate the energy separation between the ground and first excited states, in the weak tunnelling limit ( κ >> 1 ) between the two minima of the double well potential.For u < −1 , using the SHO approximation developed in Eqs. ( 9) and ( 10), we express the phase angle η in term of the eigenenergy E: where ±x 1 and ±x 2 are the turning points.In the weak tunnelling limit ( κ >> 1 ), the solutions of (20) occur near the poles of the tangent function η n ≈ (n + 1/2)π , where n = 0, 1, 2, • • • .For the ground and low-lying excited states, we can write Combining Eqs. ( 21) and ( 22), we determine the energy spacing dE 10 between the ground and first excited states ( n = 0): In the exponent, κ , the energy parameter can be approximated by that of the ground state.

Numerical results
We have discussed the analytical expressions of x c and E g using a continuum approximation, dE 10 using the WKB approximation.We now compare our continuum analytical estimates with the exact numerical solution for the spectrum of the discrete equation (2).Solid lines in Fig. 2a depict analytical results for x c .The red dots are the exact numerical results.Similar comparisons are provided for the ground state energy E 0 of the original Hamiltonian (1) and energy separation dE 10 between the ground state and first excited state (see Fig. 2b and c respectively).Clearly, the WKB estimates are very close to the exact numerical results.The exact numerical solution for dE 10 also indicates that dE 10 decays with U/t much faster for large N than for small N (Fig. 2d).This has important consequences for the decoherence and instability of "Cat states" consisting of a large number of particles, with the addition of extremely small thermal fluctuations.Figure 3a shows the energy behaviour with U at N = 50 .A double well potential emerges if U/t < −0.04 , its critical value.The separation, dE 10 , between the ground state energy ( E 0 ) and the first excited energy ( E 1 ) decreases exponentially with |U/t|, at fixed N, for u < −1.0 (see Fig. 3b).The corresponding exponential decrease of dE 10 with N, at fixed |U/t| is apparent in Fig. 2c.It is this behaviour that leads to exponentially rapid decoherence and fragility of Schrödinger Cat states.We plot (see Fig. 4) the square of the wave function, both above and below the critical point u = − 1.0 , for N = 50 particles.For N even, it is convenient to introduce the translated discrete index l t ≡ l − N/2 .Figure 4a shows the square of ground state wave function and first excited wave function at u = −0.9 .(� 1 l t ) 2 is zero at l t = 0 and has two separate maximum around the origin, whereas (� 0 l t ) 2 has one maximum at l t = 0 .For u < −1.0 , there are two Gaussian peaks in the ground state (see Fig. 4b), centered near the potential energy minima.The wings of these two Gaussian peaks superimpose around l t = 0 .Consequently, (� 0 l t ) 2 at l t = 0 is nonzero (see figure 4b).The first excited state of a symmetric potential is antisymmetric in nature.As a result (� 1 l t ) 2 is zero at l t = 0 .Figure 4b shows that the overlap between (� 0 l t ) 2 and (� 1 l t ) 2 at u = −1.25 .The separation between the two peaks becomes more pronounced with decreasing u (see Fig. 4c).Moreover, the probability densities (� 0 l t ) 2 and (� 1 l t ) 2 become almost indistinguishable from each other.This has serious implication for the off-diagonal elements of the system's thermal density matrix, which provides a measure of quantum coherence between the two components of the Schrödinger Cat state.For choices of u well below the critical point ( u = − 1.0 ), small thermal fluctuations will mix the nearly degenerate ground and first excited state.Given the close resemblance of (� 0 l t ) 2 and (� 1 l t ) 2 , there is nearly complete phase cancellation of the off-diagonal coherence in the density matrix at temperature scales corresponding to the nearly vanishing energy scale dE 10 (see Fig. 3b).

Region for robust Cat state
As described above, the ground state has two Gaussian peaks for u < − 1.0 .The separation (dl) between the two Gaussian peaks depends on U/t and N. Figure 5 shows the behaviour of dl/N in the U/t − N plane.The upper dotted line represents the critical condition u = − 1.0 for the onset of the double well potential.The lower dashed line demarks the region where a well developed "Cat state" appears.This is defined by the condition that Here, the blue dots represent the exact numerical solution of Eq. ( 2) and the solid red line is the approximate WKB analytical solution given in Eq. (7).Clearly, the WKB approximation provides an excellent estimate.For U/t = − 0.05 , a double well potential occurs if N > 40 .In panel (b), we compare the precise numerical results (blue dots) of Eq. ( 2) with the approximate WKB solutions (solid red lines) for the ground state energy E 0 of the original eigenvalue equations ( 2 |� 0 (l peak )| 2 ≤ 0.1 , where ±l peak represent the locations of the Gaussian peaks in the double-well potential.Above this line, the Gaussian probability densities overlap noticeably.Below this line, the probability densities associated with the two components of the Cat are well-distinguished.When dl/N → 1 , the ground state approaches an extreme Cat state (NOON state) of the form (|N, 0� + |0, N�)/ √ 2.

Density matrix and decoherence temperature
We now consider the finite-temperature density matrix for fixed N. The off-diagonal coherence peaks of this matrix provide a measure of the stability of Cat states to thermal fluctuations.The density matrix of this system is given by: ρ = N i=0 P(E i )|� (i) ��� (i) | , where P(E i ) = 1 Z N exp(−βE i ) is the probability distribution, β = 1 k B T , T is the temperature and Z N is the partition function in canonical ensemble.For N bosons, the Hilbert space has N + 1 basis states.Therefore, the density matrix is an (N + 1) × (N + 1) matrix.The l, l ′ th element of the matrix is given by: ρ Here both l and l ′ range from 0 to N. At low temperatures ( k B T < dE 10 ), the density matrix exhibits four peaks-two diagonal and two off-diagonal, dominated the ground state probability P(E 0 ) .As the temperature increases, more excited states contribute to the density matrix.The resulting phase cancellation from the excited states leads to a decrease in the amplitudes of the off-diagonal peaks with temperature.The ground state is a symmetric wave function about the point l = N/2 in the |l� ≡ |l, N − l� basis.At zero temperature, the density matrix has four equal peaks at four locations.The diagonal peaks in ρ l,l ′ describe "populations", whereas the off -diagonal peaks describe "coherences".At finite temperatures, the density matrix acquires contributions from the excited states, most notably the antisymmetric first excited state.For N even, it is convenient to introduce the translated integer index l t ≡ l − N/2 .Antisymmetry implies that the off- diagonal coherence products � (1) −l t are negative.These contribute to phase cancellation of positive ground state terms � (0) As a result, the thermal admixture of the antisymmetric first excited state contributes significantly to decoherence of the Cat states.We consider all exact excited states together with the ground state in our numerical calculations.
In the off-diagonal region of the ( l t , l ′ t ) grid the antisymmetric excited state coherences cancel those of the symmetric states.Consequently, the amplitude of the off-diagonal peaks in the density matrix decrease rapidly www.nature.com/scientificreports/with increasing temperature as shown in Fig. 6, panels a-d.We define the decoherence temperature ( T d,1/2 ) at which the amplitude of the off-diagonal peak is 1/2 times the amplitude of the diagonal peak at zero temperature.More generally, we define T d,R according to the condition that The ground state has two Gaussian peaks only if u < −1 .Below this critical value, there are four peaks in the density matrix at zero temperature.However, the peaks are not well separated close to u = − 1.0 .We refer to this quantum state for which the Cat components are not well-developed as a weak Cat state.Panels e-h of figure (6) illustrate the behavior of the four peaks of the density matrix with increasing temperature for the weak Cat.In this case, it is very difficult to distinguish between the peaks.
We now focus on Cat states with well-separated components, and R cat ≤ 0.1 .Figure 7 illustrates the rapid lowering of the decoherence temperature as Cat states become better-developed and Cat components are wellseparated.The color bars in Fig. 7   Table 1.Decoherence temperatures at t = 0.1 eV ( ≡ 1162 K) for N = 25 bosons.Clearly, the decoherence temperatures closely track the energy separation dE 10 between the ground and first excited states.Increasing u values above − 1.28 will increase both R cat and T d,R values.However, a large R cat > 0.5 indicates that the peaks of the density matrix and the components of the Cat state are not distinctly different.1 and 2 with the choice of t = 0.1 eV for N = 25 and N = 50 , respectively.Clearly, decoherence occurs at lower temperature for larger values of |U/t| where the Cat states are more well-defined.For a given value of u < −1.0 , the loss of coherence is even more severe for large N than for large |U/t|.The dynamical time scales describing decoherence and thermalization may be of interest in the future research 48 .

Particle number fluctuations
Another potential source of decoherence of the Schrödinger Cat state is through exchange of particles with a reservoir.We now discuss the two-site Bose-Hubbard model in the grand canonical ensemble.Here, the average particle number is determined by a chemical potential.The canonical partition function and probability for a particular N are given by: ) and P N = Z N /Z , respectively.The grand partition function is Z = N Z N .
The attractive Bose-Hubbard model ( 1) is physically unstable when placed in contact with a reservoir from which particles can enter the system.With a purely attractive, on-site Hubbard interaction, it is energetically favourable for an arbitrarily large number of particles to enter the system from the reservoir.This unmitigated accumulation of the particles is unphysical.In order the regulate this divergence, we introduce an additional repulsive term in the original Hamiltonian that is active when the particle number becomes very large.When N is not too large, the onsite attractive interaction remains dominant.The additional repulsive interaction is of the form α(n a + n b ) γ , where α is a positive energy parameter.Equation (2) in the |l� = |l, N − l� basis, then becomes: We numerically solve the above equation to determine all energy eigenvalues and eigenfunctions, for a large variety of choices of N. Figure 8a shows the behaviour of the probability distribution for different values of α at γ = 4 .Clearly, P N is zero in the large N limit and unlimited accumulation of particles from the reservoir is prevented.The value of N for which P N reaches its maximum is denoted as N max .N max depends on all other variables (α, γ , U, t, µ) .We plot P N , as function of the chemical potential µ , in µ − N plane (see Fig. 8b).The sharp jump of the peak position ( N max ) for specific µ values is straightforward to interpret.P N depends on the behaviour of the individual energy eigenvalues, E i with N. Figure 9 shows the behaviour of the ground and first excited state energies for different values of α and µ .For α = 0 , E ′ i (≡ E i − µN) first increases with N and then decreases without bound (see Fig. 9a).As discussed, above, this situation is unphysical.With addition of the term αN γ ( α > 0, γ = 4 ) to the Hamiltonian, the energy eigenvalues first increase and then possibly decrease with N, but eventually increase without bound for sufficiently large N as shown in Fig. 9b and c. Typically, E ′ i has a local minimum near N = 0 , but a second deeper local minimum appears as a function of N for µ larger than a critical value (Fig. 9c).Accordingly, we see a sharp jump in N max as µ increases.The location of the jump depends on the values of α and γ .For example, when α = 4 × 10 −6 and γ = 4 , the critical chemical potential is µ = −1.25.
The above analysis (see Fig. 8b) provides a convenient upper cutoff in particle number in our numerical simulations of the grand partition function for a given parameter range.For example, if α = 4 × 10 −6 and γ = 4 the upper cutoff of particle number can be chosen as 70.This upper cutoff represents the minimum value of N in the grand partition sum above which our numerical results are independent of the truncation of the sum over N. We use this simplification to calculate the average number of particles ( < N >= 1 Z N NZ N ) and the root mean square of particle fluctuation ( δN = √ < N 2 > − < N > 2 ) in the grand canonical ensemble.Figure 10a shows the behaviour of average number of bosons as a function of µ at U/t = −0.05 ,for different values of α .< N > is zero for large negative values of µ .As µ increases particles are eventually admitted into the system from the reservoir with a dramatic increase in < N > at a critical µ value.< N > continues to increase with further increase in µ .The values of < N > increase with the decrease of α .For example, at µ = 0 and U/t = −0.05 ,< N >= 36 for α = 10 −5 and < N >= 55 for α = 4 × 10 −6 .For the same U/t, the effective potential is a double well with two minima for N > 40 .The corresponding density matrix in the grand canonical ensemble exhibits four peaks when < N > exceeds 40.For µ < 0 and α = 10 −5 , Cat states do not exist.For α = 4 × 10 −6 , Cat states exist for −0.95 ≤ µ ≤ 0 .The root mean square fluctuation δN is shown in Fig. 10b and δN <N> is depicted in Fig. 10c.As expected, δN exhibits a sharp peak near the critical value of µ .The density matrix in the grand    Here, T g d,1/2 and T c d,1/2 are the decoherence temperatures in grand canonical and canonical ensembles, respectively.These decoherence temperatures correspond to a factor of 2 decrease in the height of the offdiagonal peaks of the density matrix.The comparison to the canonical ensemble (final column) is made by choosing the particle number N equal to the average value < N > in the grand canonical ensemble.The parameters used in the simulation are U/t = −0.05,α/t = 4 × 10 −6 , and γ = 4. www.nature.com/scientificreports/canonical ensemble is defined as ρg (l, l ′ ) = N P N ρc N (l, l ′ ) , where ρc N (l, l ′ ) is the N−particle density matrix in canonical ensemble.Here we set ρc N (l, l ′ ) = 0 for any l, l ′ > N .We calculate the decoherence temperatures in the grand canonical ensemble and compare them with those in the canonical ensemble.Table 3 lists the decoherence temperatures T c d,1/2 and T g d,1/2 in the canonical and grand canonical ensembles, respectively.The table shows decoherence temperatures in the canonical ensemble are slightly higher than those in the grand canonical ensemble.As expected, particle number fluctuations lead to further decoherence of the Cat states.However, thermal fluctuations are considerably more destructive to the coherence of the Schrödinger Cat.

Discussion and conclusions
In summary, we have shown that Schrödinger Cat states of the two-site Bose-Hubbard model are highly susceptible to loss of quantum coherence when placed in contact with a physical environment.In general, the more distinguished the components of the Cat state, the more unstable it is to decoherence.In the model we studied, this arises from the near degeneracy of the symmetric and antisymmetric quantum superpositions of the Cat components.This illustrates the elusive nature of stable quantum entanglement involving large numbers of material particles.
One possible realization of Schrödinger Cat states is in ultracold atomic systems in an external magnetic field, in which the Feshbach resonances can be utilized to control the attractive interaction between atoms 20,49 .A bosonic atom may be trapped spatially in nearby minima of a shallow optical trap that allows tunnelling between the local minima.This would require nano-Kelvin temperature scales.Another interesting possibility is in semiconductor quantum wells containing stable, bound electron-hole pairs.These bosonic excitations can couple strongly to optical cavity modes forming exciton-polaritons that have effective masses that are 10 −4 to 10 −5 times the bare electron mass.In the context of photonic crystal cavity modes, it has been suggested that these exciton-polaritons could exhibit Bose-Einstein condensation at room temperature 50,51 .When the relevant semiconductor quantum well is sandwiched by 3D photonic band gap materials, above and below, there are degenerate valleys in momentum space where the exciton polaritons may condense.In exciton systems, Feshbach resonances can occur without recourse to an external magnetic field 24 , enabling attractive interactions between the exciton-polaritons.Hopping between the degenerate valleys in the momentum space may be facilitated by interactions with phonons in the quantum well.
The realization of stable quantum superposition states of material particles on a mesoscopic scale would be a significant advance in quantum science.We hope that our analysis provides useful insights into the fundamental challenges involved and the reasons why such quantum superpositions are particularly delicate.

Figure 1 .
Figure 1.Panels (a) and (b) display the characteristics of the effective potential and the ground state wave function, respectively, for two different values of u, one above ( u = − 0.5 or U = − 0.01 ) and one below ( u = −1.25 or U = −0.025 ) the u critical value (u = − 1) .HereV 0 = UN 2 ( N 2 − 1).The other parameters used for this calculation are N = 100 and t = 1.0 .For u = − 0.5 , the ground state wavefunction has a single Gaussian peak.For u < − 1.0 , the effective potential has two local minima and the ground state wavefunction exhibits two Gaussian peaks.Panel (c) shows the classical turning points ( x 1 and x 2 ) for E = E g (the ground state energy) and the location, x c , of the potential minimum.Here, the green line and blue line denote the ground state energy and effective potential, respectively.

Figure 2 .
Figure 2. Panel (a) shows the behavior of the location of the potential minimum, x c , with N at U/t = − 0.05 .Here, the blue dots represent the exact numerical solution of Eq. (2) and the solid red line is the approximate WKB analytical solution given in Eq. (7).Clearly, the WKB approximation provides an excellent estimate.For U/t = − 0.05 , a double well potential occurs if N > 40 .In panel (b), we compare the precise numerical results (blue dots) of Eq. (2) with the approximate WKB solutions (solid red lines) for the ground state energy E 0 of the original eigenvalue equations (2) and (3).Panel (c) depicts the energy difference between the ground and first excited states dE 10 provided by Eqs.(11) and(23).Panel (d) shows behaviour of dE 10 with U/t for different values of N. The decay of dE 10 with U/t is much faster for large N than for small N. Here, all points are obtained from exact numerical solution of Eq. (2).
Figure 2. Panel (a) shows the behavior of the location of the potential minimum, x c , with N at U/t = − 0.05 .Here, the blue dots represent the exact numerical solution of Eq. (2) and the solid red line is the approximate WKB analytical solution given in Eq. (7).Clearly, the WKB approximation provides an excellent estimate.For U/t = − 0.05 , a double well potential occurs if N > 40 .In panel (b), we compare the precise numerical results (blue dots) of Eq. (2) with the approximate WKB solutions (solid red lines) for the ground state energy E 0 of the original eigenvalue equations (2) and (3).Panel (c) depicts the energy difference between the ground and first excited states dE 10 provided by Eqs.(11) and(23).Panel (d) shows behaviour of dE 10 with U/t for different values of N. The decay of dE 10 with U/t is much faster for large N than for small N. Here, all points are obtained from exact numerical solution of Eq. (2).

Figure 3 .
Figure 3. Panel (a) shows the behavior of several low energy states with U at N = 50 .In panel (b), the separation ( dE 10 ) between the ground state and the first excited state is discernible for U > 0 , but drops precipitously to zero for U < 0 .For N = 50 , the critical value of U/t is −0.04 as indicated by A in panel 3b.We label another two points one above (B) and the other below (C) the critical point u = − 1.0.

Figure 4 .
Figure 4. Panels (a-c) depict the square of the ground state and first excited state wave functions for three different values of u. (� 0 ) 2 and (� 1 ) 2 are markedly distinct for u > −1 , the critical value (see (a)).Figure (b) shows the numerical results for u = − 1.25 .Now both (� 0 ) 2 and (� 1 ) 2 have a minimum point at l t (≡ l − N/2) = 0 , but are discernibly distinct.A further decrease in u to − 1.5 causes distinction to become indiscernible (see (c)).This leads to nearly vanishing off-diagonal coherence in the density matrix at exponentially small temperature scales.
Figure 4. Panels (a-c) depict the square of the ground state and first excited state wave functions for three different values of u. (� 0 ) 2 and (� 1 ) 2 are markedly distinct for u > −1 , the critical value (see (a)).Figure (b) shows the numerical results for u = − 1.25 .Now both (� 0 ) 2 and (� 1 ) 2 have a minimum point at l t (≡ l − N/2) = 0 , but are discernibly distinct.A further decrease in u to − 1.5 causes distinction to become indiscernible (see (c)).This leads to nearly vanishing off-diagonal coherence in the density matrix at exponentially small temperature scales.

Figure 5 .
Figure 5. Depicted is the separation (dl) between the two Gaussian peaks (components of the Cat) in the ground state wave function 0 .We divided l by N to fit all results on the same scale.l/N ranges from 0 to +1 .The upper dotted line represents the locus of the critical points for which u = − 1.0 .The lower dotted line represents the locus of points below which our condition for a well-developed Cat state, R cat = 0.1 , is satisfied (see main text).

Figure 6 .
Figure 6.Panels (a-d) depict the density matrices of well-developed ( R cat = 0.09 ) Cat states at four distinct temperatures for U/t = − 0.05, N = 50 ( u = − 1.25 ).In panel (a), the density matrix contains four nearly equal-amplitude peaks at very low temperatures T << T d,1/2 .Panels (b) and (c) depict density matrices for T = T d,1/2 and T = T d,1/10 , respectively.Here, T d,R is defined as the temperature where the ratio between the diagonal peak to the off-diagonal peak is R.The height of the off-diagonal peaks becomes smaller as the temperature increases.In panel (d), T >> T d,1/10 and the amplitude of the diagonal peaks is noticeably diminished as the peaks broaden and merge.Here, the Cat components are less separated and coherence is lost.Panels (e-h) exhibit the density matrix for four distinct temperatures at U/t = − 0.044, N = 50 ( u = − 1.11 ).Since u is close to the critical point ( u = − 1.0 ), the Cat components are not well separated ( R cat = 0.98 ) and the states are referred to as weak Cat states.
depict the decoherence temperatures T d,1/2 and T d,1/10 in U/t − N plane.The blue dashed line represents the locus of points for which |� 0 (l t =0)| 2 |� 0 (l peak )| 2 = R cat , takes on the value 0.1 in 7a and 7b.The decoherence temperature is always highest near blue dashed line but drops exponentially below it.Below the blue dashed line, Cat states are more well-defined with better separated components, but the separation dE 10 between the ground and first excited state becomes exponentially smaller.As a result the Cat states lose their coherence at exponentially low temperature.

Figure 7 .
Figure 7. Decoherence temperatures ( T d,R ) for Schrödinger Cat states in N − U/t plane.Here T d,R is defined by the condition ρ(T d,R ) ρ(T=0) | off-diagonal-peak = R .In panel (a), R = 1/2 whereas in panel (b), R = 1/10 .The decoherence temperatures drop very rapidly as we move toward more well-developed Cat states.In panels (a) and (b), the blue dashed depicts the locus of points for which R cat = 0.1 , where R cat ≡ |� 0 (l t =0)| 2 |� 0 (l peak )| 2 .Decoherence temperatures are represented by the colour bar in log [k B T d,R /t] scale, where t is the hopping matrix elements between sites a and b.It is evident that Cat states are more robust and stable to small thermal fluctuations for large |U/t| and small N.This can be compared to the behaviour of dE 10 in Figs.2c and 2d.For the same values of |UN/(2t)|, Cat states with large N and small |U/t| exhibit lower decoherence temperatures.

Figure 8 .
Figure 8. Panel (a) shows the particle number probability distribution for three different values of α at γ = 4 .P N is almost zero at large N for α > 0 .The parameters used in our exact diagonalization of the many-body Hamiltonian and calculation of P N are U/t = −0.05,γ = 4, k B T/t = 0.1 and µ = −0.1.Panel (b) shows the distribution of P N in the µ − N plane.Here, P N is practically zero N > 70 .As a result, the grand partition sum is accurately described by keeping terms only up to N = 70 .The parameters used in panel (b) are U/t = − -0.05, γ = 4, α/t = 4 × 10 −6 and k B T/t = 0.01.

Figure 9 .Table 3 .
Figure 9. Ground and first excited state energies as a function of particle number for different choices of the repulsive energy coefficient α .Panel (a) shows the (unphysical) instability of the ground state and first excited state energies in the absence of regulation ( α = 0 ).Here energy decreases without bound for large N.For α > 0 , (panels (b) and (c)), the probability, P N , of N particles being admitted to the system from the reservoir tends to zero for large N.The energy eigenvalues have a global minimum at N = 0 (panel (b)) for small µ .But the global minimum shifts from N = 0 to finite N above a critical value of µ (panel (c)).As a result, there is a sharp jump in P N at the critical chemical potential.The parameters used for this simulations are (a) U = −0.05,t = 1, µ = −1.4,α = 0 , (b) U = −0.05,t = 1, γ = 4, α = 4 × 10 −6 , µ = −1.4 and (c) U = −0.05,t = 1, γ = 4, α = 4 × 10 −6 , µ = −1.2.

Figure 10 .
Figure 10.Panel (a) shows the sharp jump in the average boson number < N > with µ in the grand canonical ensemble for three different values of the repulsion energy parameter α in the attractive Bose-Hubbard model with U/t = −0.05 .For α = 4 × 10 −6 , the system admits particles from the reservoir above the critical value of µ = −1.25 .The values of < N > increase with the decrease of α .In the canonical ensemble, Cat states are possible at T = 0 for N > 40 .In the grand canonical ensemble < N >< 40 when α = 10 −5 for the range of µ shown.For α = 4 × 10 −6 , < N >= 55 at µ = 0 .In this case a Cat type solution exist in the range −0.95 ≤ µ ≤ 0 .The root mean square particle fluctuation δN and δN <N> are depicted in panels (b) and (c) respectively.The parameters used in this simulations are U/t = −0.05,γ = 4 and k B T/t = 0.1.